function FMatrix=mfccf(num,s,Fs)
%Fecha 10.10.08
%Syntax:        M=mfccf(num,s, Fs);

%Calcula y devuelve los coeficientes de una señal de voz
%s:     es la señal de entrada
%num:   es el número de coeficientes requeridos. 
%Fs:    frecuencia de muestreo
%Utiliza la función 'melbankm' del toolbox Voicebox de Mike Brooks copyright(c)
%1997 (GNU General Public License), disponible en internet libre y gratuito 
%para implemetar el banco triangular de filtros de mel (triangular mel filter bank)
%Modificaciones: Paula Sanz Leon


n   =   512;              %Número de puntos de FFT 
Tf  =   0.030;            %Duración del frame en segundos
N   =   Fs*Tf;            %Número de muestras por frames
fn  =   20;               %Número de filtros mel o ventanas
l   =   length(s);        %Número total de muestras en la señal de entrada
Ts  =   0.01;             %Paso del frame en segundos
FrameStep   =   Fs*Ts;    %Paso del Frame en muestras 
a   = 1;
b   = [1, -0.97];         %a y b son coeficientes de filtro pasa alto
noFrames    =   floor(l/FrameStep);                     %Máximo número de frames en la muestra de voz
FMatrix     =   zeros(noFrames-2, num);                 %Matriz para guardar los coeficientes
lifter      =   1:num;                                  %Elevador del índice del vector
lifter      =   1+floor((num)/2)*(sin(lifter*pi/num));   %Versión de elevador de seno aumentado. Peso de los coef.

if mean(abs(s)) > 0.01
        s=s/max(s);       %Normalizacion de amplitud.
end


%Segmentacion de la señal en frames sobrepuestos y calcula los MFCC

for	i=1:noFrames-2;
    
        frame   =s((i-1)*FrameStep+1:(i-1)*FrameStep+N);    %Mantiene los frames individuales. Frame actual.
        Ce1     =sum(frame.^2);                             %Energía del frame
        Ce2     =max(Ce1,2e-22);                            %floors to 2 X 10 raised to power -22
        Ce      =log(Ce2);
        framef  =filter(b,a,frame);                         %Filtro pasa alto pre-enfasis
        F       =framef.*hamming(N);                        %multiplica cada frame con una ventana hamming
        plot(frame)
        FFTo    =fft(F,N);                                  %calcula la fft del frame
        melf    =melbankm(fn,n,Fs);                         %crea n filtros, mel filter bank para el frame
        halfn   =1+floor(n/2);    
        spectr1 =log10(melf*abs(FFTo(1:halfn)).^2);         %el resultado es filtrado en la escala de mel
        spectr  =max(spectr1(:),1e-22);
        c       =dct(spectr);                               %obtiene DCT,Transformada de Coseno Discreta
                                                            %cambia dominio cepstral
        c(1)            =Ce;                                %reemplaza el primer coeficiente
        coeffs          =c(1:num);                          %retiene los primeros num coeficientes
        ncoeffs         =coeffs.*lifter';                   %
        FMatrix(i, :)   =ncoeffs';                          %asigna los mfcc a sucesivas filas i
end

%Llama a la función deltacoeff para calcular las derivadas de los coeficientes MFCC
%Concatena las matrices de derivadas para dar una matriz con 3*num columnas 
% 
% d       =(deltacoeff(FMatrix)).*0.6;    %Calcula derivada de mfcc
% d1      =(deltacoeff(d)).*0.4;          %segunda derivada
% FMatrix =[FMatrix,d,d1];                %concatena todo junto
